import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy.optimize import minimize
import os
from datetime import datetime, timedelta

# 设置中文显示
plt.rcParams['font.sans-serif'] = ['SimHei']
plt.rcParams['axes.unicode_minus'] = False

EXCEL_PATH = "E:/M/MathModel/23国赛C题/ds/历年6月底到7月初.xlsx"
OUTPUT_DIR = "E:/M/MathModel/23国赛C题/ds/利润优化结果"



def load_and_prepare_data(file_path):
    df_2020 = pd.read_excel(file_path, sheet_name='2020')
    df_2021 = pd.read_excel(file_path, sheet_name='2021')
    df_2022 = pd.read_excel(file_path, sheet_name='2022')
    df_2023 = pd.read_excel(file_path, sheet_name='2023')

    # 合并数据
    df = pd.concat([df_2020, df_2021, df_2022, df_2023])

    # 转换日期格式
    df['销售日期'] = pd.to_datetime(df['销售日期'])

    # 修正分类名称
    df['分类名称'] = df['分类名称'].str.replace('食用菌类', '食用菌')

    return df


def calculate_historical_parameters(df):
    # 按品类分组计算平均批发价和平均损耗率
    params = df.groupby('分类名称').agg(
        avg_wholesale_price=('批发价', 'mean'),
        avg_loss_rate=('损耗率(%)', 'mean')
    ).reset_index()

    return params


def get_regression_parameters():
    regression_params = {
        '花菜类': {
            'workday': (-0.0438, 11.369),
            'non_workday': (-0.0438, 12.177)
        },
        '花叶类': {
            'workday': (0.0187, 3.169),
            'non_workday': (0.021, 2.179)
        },
        '辣椒类': {
            'workday': (0.0142, 7.921),
            'non_workday': (0.0142, 7.921)
        },
        '茄类': {
            'workday': (-0.0221, 7.397),
            'non_workday': (-0.0221, 7.397)
        },
        '水生根茎类': {
            'workday': (0.006, 15.337),
            'non_workday': (0.006, 15.337)
        },
        '食用菌': {
            'workday': (0.0376, 10.896),
            'non_workday': (0.0376, 10.896)
        }
    }
    return regression_params


def is_workday(date):
    return date.weekday() < 5


def profit_function(x, categories, dates, regression_params, wholesale_prices, loss_rates):

    n_categories = len(categories)
    n_days = len(dates)
    total_profit = 0

    # 解析决策变量
    idx = 0
    for i, category in enumerate(categories):
        wholesale_price = wholesale_prices[category]
        loss_rate = loss_rates[category] / 100  # 转换为小数

        for j, date in enumerate(dates):
            P_ij = x[idx]  # 售价
            X_ij = x[idx + 1]  # 销售量
            T_ij = x[idx + 2]  # 补货量
            idx += 3

            # 计算变动成本
            B_ij = loss_rate * P_ij * T_ij

            # 计算每日利润
            daily_profit = P_ij * X_ij - B_ij - T_ij * wholesale_price
            total_profit += daily_profit

    # 返回负值用于最小化
    return -total_profit


def inventory_constraints(x, categories, dates):
    n_categories = len(categories)
    n_days = len(dates)
    constraints = []

    # 解析决策变量并创建约束
    idx = 0
    for i, category in enumerate(categories):
        # 初始库存为0
        inventory = 0

        for j, date in enumerate(dates):
            # 获取决策变量
            P_ij = x[idx]  # 售价
            X_ij = x[idx + 1]  # 销售量
            T_ij = x[idx + 2]  # 补货量

            # 计算可用库存 = 前日库存 + 今日补货
            available_inventory = inventory + T_ij

            # 约束: 销售量 <= 可用库存
            constraints.append(available_inventory - X_ij)

            # 更新库存 = 可用库存 - 销售量
            inventory = available_inventory - X_ij

            # 移动到下一组变量
            idx += 3

    return np.array(constraints)


def regression_constraints(x, categories, dates, regression_params):
    n_categories = len(categories)
    n_days = len(dates)
    constraints = []

    # 解析决策变量并创建约束
    idx = 0
    for i, category in enumerate(categories):
        for j, date in enumerate(dates):
            # 获取决策变量
            P_ij = x[idx]  # 售价
            X_ij = x[idx + 1]  # 销售量

            # 确定工作日/非工作日
            workday_type = 'workday' if is_workday(date) else 'non_workday'

            # 获取回归参数
            a_ij, b_ij = regression_params[category][workday_type]

            # 约束
            constraints.append(P_ij - (a_ij + b_ij * X_ij))

            # 移动到下一组变量
            idx += 3

    return np.array(constraints)


def non_negative_price_constraint(x, categories, dates, regression_params):
    n_categories = len(categories)
    n_days = len(dates)
    constraints = []

    # 解析决策变量并创建约束
    idx = 0
    for i, category in enumerate(categories):
        for j, date in enumerate(dates):
            # 获取决策变量
            X_ij = x[idx + 1]  # 销售量

            # 确定工作日/非工作日
            workday_type = 'workday' if is_workday(date) else 'non_workday'

            # 获取回归参数
            a_ij, b_ij = regression_params[category][workday_type]

            # 约束
            constraints.append(a_ij + b_ij * X_ij)

            # 移动到下一组变量
            idx += 3

    return np.array(constraints)


def optimize_profit(categories, dates, regression_params, wholesale_prices, loss_rates):
    n_categories = len(categories)
    n_days = len(dates)
    n_vars = n_categories * n_days * 3  # 每天每品类3个变量(P, X, T)

    x0 = []
    for category in categories:
        for date in dates:
            workday_type = 'workday' if is_workday(date) else 'non_workday'
            a, b = regression_params[category][workday_type]

            # 合理的初始值：确保价格非负
            est_sales = max(1, -a / b) * 1.2 if a < 0 else 10
            est_price = max(0.1, a + b * est_sales)  # 确保价格为正
            est_replenishment = est_sales * 1.1  # 补货略大于销量
            x0.extend([est_price, est_sales, est_replenishment])
    x0 = np.array(x0)

    # 设置更合理的边界
    bounds = []
    for i, category in enumerate(categories):
        for j, date in enumerate(dates):
            bounds.append((0.5, 30))
            bounds.append((50, 500))
            bounds.append((50, 500))

    constraints = [
        # 库存约束
        {'type': 'ineq', 'fun': lambda x: inventory_constraints(x, categories, dates)},

        # 回归关系约束
        {'type': 'eq', 'fun': lambda x: regression_constraints(x, categories, dates, regression_params),
         'tol': 1e-3},  # 增加容差

        # 价格非负约束
        {'type': 'ineq', 'fun': lambda x: non_negative_price_constraint(x, categories, dates, regression_params)}
    ]

    options = {'disp': True, 'maxiter': 500, 'ftol': 1e-4}

    solution = minimize(
        lambda x: profit_function(x, categories, dates, regression_params, wholesale_prices, loss_rates),
        x0,
        method='SLSQP',
        bounds=bounds,
        constraints=constraints,
        options=options
    )

    if not solution.success:
        print(f"优化警告: {solution.message}")
        # 即使不完全成功，也尝试返回当前最优解
        return solution, solution.x

    return solution, solution.x


def parse_optimal_values(optimal_values, categories, dates):
    n_categories = len(categories)
    n_days = len(dates)

    # 创建结果DataFrame
    results = []
    idx = 0

    for i, category in enumerate(categories):
        for j, date in enumerate(dates):
            # 获取决策变量
            P_ij = optimal_values[idx]  # 售价
            X_ij = optimal_values[idx + 1]  # 销售量
            T_ij = optimal_values[idx + 2]  # 补货量
            idx += 3

            results.append({
                '品类': category,
                '日期': date,
                '售价(元/千克)': P_ij,
                '销售量(千克)': X_ij,
                '补货量(千克)': T_ij
            })

    return pd.DataFrame(results)


def calculate_daily_profit(df_optimal, wholesale_prices, loss_rates):
    profits = []

    for _, row in df_optimal.iterrows():
        category = row['品类']
        price = row['售价(元/千克)']
        sales = row['销售量(千克)']
        replenishment = row['补货量(千克)']

        # 获取批发价和损耗率
        wholesale_price = wholesale_prices[category]
        loss_rate = loss_rates[category] / 100  # 转换为小数

        # 计算变动成本
        B_ij = loss_rate * price * replenishment

        # 计算每日利润
        daily_profit = price * sales - B_ij - replenishment * wholesale_price
        profits.append(daily_profit)

    df_optimal['利润(元)'] = profits
    return df_optimal


def main():
    # 创建输出目录
    os.makedirs(OUTPUT_DIR, exist_ok=True)

    # 加载和预处理数据
    print("加载和预处理数据")
    try:
        df = load_and_prepare_data(EXCEL_PATH)
        print(f"加载成功，共 {len(df)} 条记录")
    except Exception as e:
        print(f"数据加载失败: {e}")
        return

    # 计算历史参数
    print("计算历史参数")
    params_df = calculate_historical_parameters(df)

    # 创建批发价和损耗率字典
    wholesale_prices = {}
    loss_rates = {}
    for _, row in params_df.iterrows():
        wholesale_prices[row['分类名称']] = row['avg_wholesale_price']
        loss_rates[row['分类名称']] = row['avg_loss_rate']

    # 获取回归参数
    print("获取回归参数")
    regression_params = get_regression_parameters()

    # 定义日期范围 (2023年7月1日-7日)
    start_date = datetime(2023, 7, 1)
    dates = [start_date + timedelta(days=i) for i in range(7)]

    # 定义品类
    categories = ['花叶类', '花菜类', '辣椒类', '茄类', '水生根茎类', '食用菌']

    # 执行利润优化
    print("执行利润优化")
    solution, optimal_values = optimize_profit(
        categories, dates, regression_params, wholesale_prices, loss_rates
    )

    if solution is None:
        print("优化失败，无法继续")
        return

    # 解析优化结果
    print("解析优化结果")
    df_optimal = parse_optimal_values(optimal_values, categories, dates)

    # 计算每日利润
    print("计算每日利润")
    df_optimal = calculate_daily_profit(df_optimal, wholesale_prices, loss_rates)

    # 保存
    print("保存结果")
    output_path = os.path.join(OUTPUT_DIR, '利润优化结果.xlsx')
    df_optimal.to_excel(output_path, index=False)
    print(f"优化结果保存至: {output_path}")

    # 可视化
    print("生成可视化图表")
    try:
        # 创建透视表用于可视化
        pivot_price = df_optimal.pivot(index='日期', columns='品类', values='售价(元/千克)')
        pivot_sales = df_optimal.pivot(index='日期', columns='品类', values='销售量(千克)')
        pivot_profit = df_optimal.pivot(index='日期', columns='品类', values='利润(元)')

        # 最优售价趋势图
        plt.figure(figsize=(14, 8))
        for category in pivot_price.columns:
            plt.plot(pivot_price.index, pivot_price[category], 'o-', label=category)

        plt.title('2023年7月1-7日最优售价')
        plt.xlabel('日期')
        plt.ylabel('最优售价 (元/千克)')
        plt.legend()
        plt.grid(True, linestyle='--', alpha=0.7)
        plt.tight_layout()

        plot_path = os.path.join(OUTPUT_DIR, '最优售价趋势图.png')
        plt.savefig(plot_path)
        plt.close()
        print(f"最优售价趋势图保存成功: {plot_path}")

        # 销售量与补货量对比图
        plt.figure(figsize=(14, 10))
        for i, category in enumerate(categories, 1):
            plt.subplot(3, 2, i)
            cat_data = df_optimal[df_optimal['品类'] == category]

            plt.plot(cat_data['日期'], cat_data['销售量(千克)'], 'o-', label='销售量')
            plt.plot(cat_data['日期'], cat_data['补货量(千克)'], 's--', label='补货量')

            plt.title(f'{category}销售量与补货量')
            plt.xlabel('日期')
            plt.ylabel('数量 (千克)')
            plt.legend()
            plt.grid(True, linestyle='--', alpha=0.7)

        plt.tight_layout()
        plot_path = os.path.join(OUTPUT_DIR, '销售量与补货量对比图.png')
        plt.savefig(plot_path)
        plt.close()
        print(f"销售量与补货量对比图保存成功: {plot_path}")

        # 每日利润分布图
        plt.figure(figsize=(14, 8))
        pivot_profit.plot(kind='bar', stacked=True, figsize=(14, 8))
        plt.title('2023年7月1-7日每日利润分布')
        plt.xlabel('日期')
        plt.ylabel('利润 (元)')
        plt.legend(title='品类')
        plt.grid(axis='y', linestyle='--', alpha=0.7)
        plt.tight_layout()

        plot_path = os.path.join(OUTPUT_DIR, '每日利润分布图.png')
        plt.savefig(plot_path)
        plt.close()
        print(f"每日利润分布图保存成功: {plot_path}")

    except Exception as e:
        print(f"图表生成失败: {e}")

    # 输出最终结果
    print("\n===== 2023年7月1-7日利润优化结果 =====")
    print(f"总利润: {-solution.fun:.2f} 元")

    # 按品类汇总结果
    category_summary = df_optimal.groupby('品类').agg(
        总销售量=('销售量(千克)', 'sum'),
        总补货量=('补货量(千克)', 'sum'),
        平均售价=('售价(元/千克)', 'mean'),
        总利润=('利润(元)', 'sum')
    ).reset_index()

    print("\n各品类汇总结果:")
    print(category_summary)

    # 按日期汇总结果
    date_summary = df_optimal.groupby('日期').agg(
        总销售量=('销售量(千克)', 'sum'),
        总补货量=('补货量(千克)', 'sum'),
        总利润=('利润(元)', 'sum')
    ).reset_index()

    print("\n每日汇总结果:")
    print(date_summary)

    print(f"\n所有结果已保存至: {os.path.abspath(OUTPUT_DIR)}")


if __name__ == "__main__":
    main()